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This paper extends previous analysis on the impact of sensing noise for the navigation 
and control aspects of formation flying spacecraft. We analyze the use of Carrier-phase 
Differential GPS (CDGPS) in relative navigation filters, with a particular focus on the filter 
correlation coefficient. This work was motivated by previous publications which suggested 
that a “good” navigation filter would have a strong correlation (i.e., coefficient near —1) 
to reduce the semimajor axis (SMA) error, and therefore, the overall fuel use. However, 
practical experience with CDGPS-based filters has shown this strong correlation seldom 
occurs (typical correlations ss —0.1), even when the estimation accuracies are very good. 

We derive an analytic estimate of the filter correlation coefficient and demonstrate that, 
for the process and sensor noises levels expected with CDGPS, the expected value will be 
very low. It is also demonstrated that this correlation can be improved by increasing the 
time step of the discrete Kalman filter, but since the balance condition is not satisfied, the 
SMA error also increases. These observations are verified with several linear simulations. 

The combination of these simulations and analysis provide new insights on the crucial role 
of the process noise in determining the semimajor axis knowledge. 

I. Introduction 

In formation flying missions, accurate knowledge of the difference in semimajor axes, or equivalently, 
the difference in orbital energy, between the vehicles in the formation is important HE[ 7 A difference in 
semimajor axes means that the two vehicles have different orbital periods and thus they will drift out of 
formation unless considerable control effort is applied. 3 The output of the CDGPS Kalman filter includes 
the relative formation state in a Local Vertical Local Horizontal (LVLH) reference frame. Understanding 
the relationship between position and velocity accuracies and semimajor axis accuracy is key to evaluating 
the output of this type of filter. While Ref. (4] develops the navigation error analysis from absolute state 
relations, the results can be reformulated for the relative case. Ref. [5] analyzed that case in detail, comparing 
linear and nonlinear simulation results with analytic predictions. This paper extends the previous work to 
consider further modifications to the Carrier-Phase Differential GPS (CDGPS) Kalman filter to determine 
their effect on the correlation and balance , which are the two key properties. This includes further details 
on the analytic predictions, an analysis of the effect of adding velocity sensing, and an investigation of the 
impact of increasing the time step in the discrete Kalman filter. 

The relative navigation error equations, shown below, relate semimajor axis error to position and velocity 
errors. Note that this discussion is limited to circular reference orbits. The standard deviation of the 
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Fig. 1 : Contours of constant semimajor axis vs. po- 
sition and velocity accuracy. Contours given for 3 
levels of correlation. 


Contours of Constant Balance for g =1 



Fig. 2 : Contours of constant balance (straight lines) 
illustrate region where correlation affects SMA accu- 
racy. Numbers indicate balance index. 


differential semimajor axis estimate, ctao., is given bj3l 


<JAa = 2\/4o-2 + - p x ycr x ay + — cr? 


(1) 


The parameters a x , <jy, and p X y are derived from the error covariance matrix for the relative LVLH state 
estimate, x = [ x y x y ] T with estimation error x = x — x, which is assumed to be unbiased, E[x] = 0 , 
and have a covariance 


E[ 


xx T l = 



Pxy® x ®y 
' 

Pxx® X ® X 

Pxy®x ®y 


Pyx&y&x 


Pyx®y ®x 
2 

PyyCJyGy 

( 2 ) 

Pxx^x O x 

Pxy® x ®y 

®x 

Pxy® x ®y 


Pyx&y&x 

PyyO’yO’y 

Pyx® y® x 

4 \ 



Note that if the radial position and in-track velocity are linearly correlated ( p X y = —1), the expression for 
semimajor axis variance from Eq. [T] reduces to 


<TAa = 2\/4a2 _ i axa . + T CT 2 = 2^ 


2a?, a. 


( 3 ) 


If the position and velocity error are linearly correlated and satisfy 


— 2na x 


( 4 ) 


then the position and velocity errors cancel and there is no error in the semimajor axis estimate. In other 
words, the two requirements for zero semimajor axis variance are: 


Pxy = — 1 and <jy = 2na x 

which will subsequently be referred to as the correlation and balance requirements. 

The relationship between a x , <7y , p X y, and UAa is illustrated in Fig. [l|4l The x and y axes of the plot 
are the standard deviations of the position and velocity estimation errors. Contours of constant semimajor 
axis standard deviation are shown on the figure. Each contour is associated with a value of p X y in addition 
to a level of OAa\ several values of p X y are shown for each level of a a a- The diagonal of peaks indicates 
where <jy = 2 na x . Along the diagonal of peaks, the lines of constant semimajor axis experience a “bump” 
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Fig. 3 : Contours of constant Q and R are shown on 
axes of position and velocity accuracy. 


Fig. 4: Contours of constant semimajor axis shown 
vs position/ velocity accuracy. Contours of constant 
Q and R are dotted lines in background. 


that increases in size as the correlation tends towards — 1. This bump corresponds to increasing cancelation 
between the error in x and y that results from increasing correlation in these errors. Essentially, if the 
errors have high correlation and the proper balance, the higher error levels can be tolerated with the same 
resulting semimajor axis error. Each point on the graph corresponds to a unique set of a x and ay. However, 
many points on the graph are intersected by more than one contour of constant semimajor axis. It is the 
correlation that determines the specific contour on which the system lies. Fig. [2] shows the contours for a 
constant value of the semimajor axis error for varying values of correlation. The plot also shows lines of 
constant balance index , which is calculated as 


bal = 


1 - 


2na x 


( 5 ) 


Clearly bal should be zero when the balance requirement is met. However, the figure shows that if the 
balance is slightly off, i.e.. bal « 0.5 or bal ss | — 1|, then the effects of the correlation are reduced. 


II. Linear Planar Model Simulations with Q and R 

This section investigates the relationship between the Kalman filter parameters and the resulting estimate 
accuracy. These results are based on the premise that the Kalman filter produces the best answer from the 
given the models and measurements, and that the balance or correlation between elements of the state vector 
estimate are not important. Initial investigations showed that, predictably, the levels of measurement and 
process noise have the most influence of the estimate accuracy. These levels are indicators of how well the 
sensors and the dynamics are modeled. The relative levels of these noises determines how the filter will weigh 
new measurement information against the current state estimate propagated from the dynamics information. 
A CDGPS navigation filter has nonlinearities in both the system dynamics and the measurement equations, 
and because the set of visible GPS satellites changes, the measurement matrix H will change, and the state 
vector length will grow or shrink as the set of estimated carrier biases changes!-^ These factors make it 
difficult to understand direct relationships between the filter parameters and the navigational accuracies. 
Thus, we started with a simplified Linear Planar Model (LPM) to develop insights into the behavior of a 
relative navigation filter using CDGPS. 5 

Meeting the balance and correlation requirements discussed in Section [I] corresponds to being on the 
“bump” in Fig. [T] but the baseline results of the LPM simulation were above this region. This leads to the 
question of how changing filter inputs will move the output closer to or further from the bump. To answer 
this question, the LPM simulation was run for a range of measurement and process noise levels. 

For each unique assignment of Q and R, the resulting error variances for radial position and in-track 
velocity, a x and ay, were recorded. The corresponding semimajor axis error, a\ a , was calculated using 
Eq. [l] Fig. [3] shows lines of constant Q and R on axes of a x and ay. The diagonal of peaks on this graph 
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indicates the location of the “bump,” where the balance and correlation requirements are met (which means 
< Ty = 2 na x ). By moving from one line of constant Q or R to another, one can see how decreasing the process 
or measurement noise would change the resulting position and velocity error. However, it is semimajor axis 
knowledge, not just position and velocity knowledge, that is important for control performance. 

Several graphs demonstrate the relationship between Q, R and o\ a . First, Fig. [4] is reproduced and 
lines of constant semimajor axis error added. The lines for constant Q and R are dimmed for clarity. Note 
the lines of constant semimajor axis are horizontal, which corresponds with the horizontal sections of the 
semimajor axis contours on Fig. |T] The effect of changes in Q and R on o\ a can be assessed by looking at the 
constant lines for all three values. Because the lines of constant semimajor axis are horizontal, improvement 
in ctao can only be accomplished by moving in the vertical direction on the graph, which is equivalent to 
decreasing <jy. Whether this requires decreasing Q or R depends on the angles between the horizontal lines 
of constant a^a and the contours of constant Q and R. 

The LPM results shown thus far occupy the region above the Oy = 2na x line. Ref. [5] investigated three 
questions about the system: (i) What features of the system cause this behavior? (ii) Does this agree with 
expected behavior? (iii) What is required to cause the LPM results to move closer to this line, and below 
the line? 

This was done by modifying the LPM to put more emphasis on the orbital dynamics in the Kalman filter. 
This was based on the observation that it is reasonable to expect the correlation between radial position and 
in-track velocity will increase only when the estimate depends more on the dynamics model embedded in 
the filter. The following extends the previous analysis to consider the effects of other changes to the Kalman 
Filter and to expand the analytic results. 

III. Discrete Simulations for Varying At 



One way to put more emphasis on the dy- 
namics model is to reduce the measurement 
rate. Prior analysis of the CDGPS filter by 
Busse 12 used a 1 Hz rate, although much longer 
time steps were considered in Ref. 11 . One dif- 
ficulty with increasing the propagation time is 
that the nonlinearity in the orbit propagation 
becomes much more significant, and more so- 
phisticated models and propagation algorithms 
must be used, especially for the error covari- 
anceliilliH The following investigation of the 
effect of varying measurement update rate by 
changing the discretization time step uses the 
linear model and thus ignores these effects. 

To do this, a family of discrete Riccati equa- 
tions was solved using Hill’s dynamics at various 
discretization times. In each case, the constant 
spectral density matrix Q c associated with the 
process noise of the continuous dynamics was 
converted to the appropriate discrete process 

noise Qd using the conversion algorithm in Ref. Il3] (DISRW in Chapter 9). The simulation results in Figure[5] 
were done for several values of a qc = { 10 — 5 , 10“°75 x lCC 7 }m/s 3 / 2 ( Q c = o\ c I- 2 ) and a constant sensing noise 
covariance R = (5 x 10~ 3 m) 2 . The plots shows that, as the discretization time step is increased, the corre- 
lation coefficient tends to — 1. This trend is expected, because longer propagation times will force increased 
filter dependence on the equations of motion, translating into increased overall correlatiorj^J However, the 
proper position to velocity error balance is not achieved, causing the SMA accuracy to degrade. The plot 
clearly shows the role of the process noise Q c in determining the SMA error growth. 


Fig. 5 : Filter correlation, balance, and SMA error for 
various discretization time steps At (seconds). 


a One caution: the correlation approaches —1 for filter time steps larger than 1000 seconds for which the nonlinearity in the 
dynamics will play an important role. 
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Comparing this result to Fig. [2] it can be seen that while the balance index is not zero, it is producing 
the effect of lowering the overall semimajor axis in combination with the correlation. However, the canceling 
effects of correlation and balance are not sufficient to prevent the semimajor axis error from growing rapidly 
as the time step is increased. Thus, decreasing the measurement update rate to produce higher correlation 
is not a viable strategy for reducing semimajor axis error. 


IV. Agreement with Analytical Work 


The analysis in this section develops analytical relationships between the process and measurement noise 
levels and the errors and correlations. The link between the two is established through the algebraic Riccati 
equation. The matrix Riccati equation has no known analytic solution for system more complicated than a 
double integrator. However, for a Kalman filter with a time step At « 1 second for a 90 minute orbit, the 
coupling between motion in the x and y directions is very weak and the dynamics can be well approximated 
as two weakly coupled double integrators. 

Start with the planar equations of motion in an LVLH frame about a circular reference orbit with 
frequency n 


x — 2 ny — 3 n 2 x = f x 
y + 2 nx = f y 


( 6 ) 

( 7 ) 


where f x and f y are disturbance accelerations with identical spectral densities <Tq . We can write this in state 
space form using X = [ x , x 1 y, y] T and dynamics matrices Ah and Bh . When only position measurements 
are available, the output matrix is 


H h = 


10 0 0 
0 0 10 


( 8 ) 


with sensing noise given by <j 2 r on each measurement. These dynamics can be transformed to a new 
system of equations with the state X = [x, x ' , y , y '] T , where (*)' = (*)Af, X = T p X, and T p = 
diag([ 1 At 1 At ]), to yield 


x" = 2 ey' + 3e 2 x + (At) 2 f x (9) 

y" = -2ex' + (At) 2 f y (10) 


with e = nAt. Note that when e <C 1, which is true for this application, the x and y dynamics can essentially 
be written as double integrators for which the solution of the Riccati equation is easily found. The combined 
dynamics then consist of these two double integrators with coupling terms of order e and e 2 . 


A. Correlation with Position Measurements 


The dynamics in Eqs. [9] and 10 are written in state space form as 

X' = AX + Bf (11) 

where A = At^TpAhTp 1 ), B = At(T p Bh), and H = Hh . The differential filter Riccati equation for the 
original system (Eqs. [6] and [7]) must also be transformed using T p , and the result is that, at steady state 

0 = AP + PA t + BQB t - PH t R~ 1 HP (12) 

where Q = (<7g/Af)/ 2 and R = (a R / At)I 2 , with the factors of 1/ At resulting from the transformation. To 
proceed, the covariance for the transformed state is represented as a Taylor expansion in e. Substituting the 
expressions for P = P 0 + tP\ + e 2 P 2 and A = A 0 + eAi + e 2 A 2 in the Riccati equation and grouping terms 
in the same power of e, it is possible to solve for the expansion of the covariance matrix (Pkxxi Pky'y', and 
Pkxy 1 j A: = 0,1, 2). The A matrix in Eq. 11 gives the expansion matrices Aq and A\ 


An — 


Aio — 


^4oo O 2 
O 2 ^4oo 

0 1 
0 0 


A, = 


A 1 1 — 


0 2 An 

— An O2 

0 0 
0 2 


(13) 


(14) 
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Define x = [ x x' ] T and y = [ y y' } T ■ Then P 0xy = 0 for the solution of the two independent double 
integrators. The first-order solution to the Riccati equation is computed using 


A = 

and is found to be 


0 1 
0 0 


B = 


0 

At 2 


H=[ 1 0 ] 


Pqxx — Pnvv — 


Oyy 




1 3 

Q a R 


aQa R At 

a Q a R At VZcr^cr^At 2 
Substituting into the Riccati equation and isolating 0(e) terms yields 

PoAi + PiAq + AqPi + A\Pq = PqDPi + P\DPq 
Pi{Al - DP 0 ) + (A 0 - P 0 D)P 1 = ~(P 0 A{ + A X P b) 

where 


D = 


c 

0 2 ’ 

g = P 

1 

0 

. °2 

c 

a R 

0 

0 


Using knowledge that Po xx = Po yy , 


where 


{K - DP 0 ) = 


A tP 0 . 


E O 2 
O 2 E 

AtP 0 , 


and 


Defining G as 


E = 


(A, - PoD) = 


G — AiPo xx - P 0xx A n — 


E t 

0 2 


0 2 

E t 


0 2 

2Tbxx' 


-2 F 0XX ' 

0 2 


(15) 

(16) 

(17) 

(18) 

(19) 

( 20 ) 
(21) 
(22) 
(23) 


(Po^T + A-Po) — 


0 2 G 
-G 0 2 


Substituting into the Riccati equation yields 

Pi x 


P\ xx 

Pit 


Pi. 


E 0 2 
0 2 E 


E t 

0 2 


0 2 

E t 


P, xx 

Pit 


which is a system of three equations 


Pi E + E t P 1 

-■-xx a x 


=0 4 


Pi ; 

Pi. 


Pi E + E 1 Pi =0 4 

1 yy 1 yy ^ 

Pi P + P T Pi = — G 

-*-xv -*-xv 


The solution to Eqs. 26 and 27 is Pi xx = Pi yy =0, and Eq. 28 gives 


Pixy' — 


2eriP 


Qxx’ 


AtPc 


Plx'y = “Pi 


1 

1 

0 

to 

1 

1 

(M 

O 

O 

l 


Oxx 


1 xy' 


Plx'y' — Plxv — 0 


lxy 


(24) 


(25) 


(26) 

(27) 

(28) 


( 29 ) 
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Then the correlation coefficient can be found from 


Prr 


Pxy — 


Pxy'/At 


eP 


1 xy' 


y/ Pxx Pyy \J PxxPy'y' / yj PoxxP()y'y' 


(30) 


Using the expressions given above, it follows that 


Pxy 


2(nAfK(, ?? .At) 

,/oA^2„2 V 


\/2/\t(jQ(J R 


a Q a R At ‘ 


- 1/2 



/ — T 3 i 1 

V%cr R 2 <Jq 2 At 

\J 2 a Q a R At ' 2 


(31) 

(32) 


which, upon substitution in Eq. [T] gives a semimajor axis error of 

2i a i 

CTAa = — CTqCTh (33) 

Eqs. [32] and [33] predict the correlation and semimajor axis error variance that can be expected from a contin- 
uous Kalman filter. Since these predictions were developed from a model based on two double integrators, 
the results may not immediately map to the full CDGPS problem. However, these equations agree well with 
LPM simulations!® 

To verify the ranges of <jq and <r R for which the equation is valid. The correlation coefficient, p X y , has a 
limit of —1, so 

1 =*> > n 2 (34) 

V a Q a R 

In low Earth orbit (n ss 0.001), Eq. |34| implies that ctr must be no more than six orders of magnitude larger 
than aQ. Many space-rated GPS receivers produce differential carrier phases measurements with millimeters 
of error. To achieve p X y = —1, the dynamics environment of the vehicle would have to be modeled to 
nanometer-level accuracy, which is not currently possible. For a typical example with ug = 1 x 10 -6 m/s 3 / 2 
and ctr = 5 x 10 -3 ms 1 / 2 , note that 


/ 5 x 10 -3 

\p X y\ « 1 X 10-y r - TF? = 0.0707 « 1 

This interpretation of the relationship between correlation and process and measurement noises concurs with 
conclusions drawn from the LPM examples. 


B. Examination of Balance Requirement 

One strategy for minimizing the semimajor axis error is to achieve a high correlation and to simultaneously 
have a balance of errors given by 

ay/a x = 2n (35) 

For the transformed state described by Eqs. [9|and[T0| this requirement is 

a y i/a x = 2nAt (36) 

Section |A| showed that for Hill’s equations, to first order, the standard deviations for in-track velocity and 
radial position estimates are 


<V « = (V2a% 2 a){ 2 At 2 ) 1/2 

*x = s/P^ « VPoZ = (\/ 2 ^ / 2 4 /2 ) 1/2 

Substituting into Eq. [36] gives the analytic expression for the balance condition 


V 


/&x = 


V^oy/V^Ai 2 


\ 


J Q U R 
1/2 3/2 

a Q VR 


= AtJ^- 
V a R 


(37) 

(38) 


( 39 ) 
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Using Eq. [36] perfect balance requires that 


A tJ— = 2nAt 


1 

2 n 


(40) 


This expression can be substituted into Eq. [32] to find the correlation that is achieved when the balance is 
correct 


I&R 


Pxy = -nj— = ~ n ( 2^ ) = _0 - 5 


a Q 


(41) 


This analysis suggests that achieving the required balance in this example is incompatible with achieving a 
correlation of — 1. This indicates that the strategy of minimizing semimajor axis by correlation and balance 
is impossible. 


C. Correlation with Position and Velocity Measurements 

In this example, the system will be augmented with a velocity sensor. The sensor noise and measurement 
matrices for the transformed double integrator dynamics now take the form 


R = 


0 


R = R/At 


H = In 


H = 


1 0 
0 At,- 1 


(42) 


and the process noise term Q = Q/ At. Substituting into the Riccati equation and solving for the terms in 
the covariance matrix yields 


Pn = V2a% 2 a)J 2 


P 22 = y/Zall 2 crl! 2 At 2 


22 
Pl2 = 


rx u Q 
G rx & 


(1 + 0.5/V) 1/2 

(1 + K) 

(1 + 0.5A^ 1/2 

(1 + 1<) 


(43) 


(1 + /0 

where K = 


Substituting Eq. [43] into the Hill’s expansion used in Section [A] yields a matrix equation in the same form 
as Eq. |25[ however, the matrix D is now redefined as 


D = 


C 0 2 

o 2 c 


c = 


At I a 2 


0 


0 l/(o-r± Ai ) 


(44) 


This set of matrix equations once again indicates that Pl^ = P\ = 0 2 . Solving the system of matrix 
equations yields 


Pi = Pi , , = 0 2 
Pi , = -Pi , = —2 


. a 2 P 0 _ 


AtPn. 


1 


Solving for the correlation coefficient, p X y , gives 

ePi„, 

Pxy 


(P^Po^)* 



1 + K 


1 

1 + 0.5AT 


(45) 

(46) 

(47) 


From Eq. 47 it can be seen that as the velocity measurement becomes more accurate ( K increases), the 
correlation magnitude is reduced. This is consistent with the observation that the filter makes more use 
of dynamics when there are fewer measurements. Substituting this value for p X y into the expression for 
semimajor axis standard deviation yields 


24 i 3 

(7 Aa = GtxOq 

n W 


'vQVrx + 4 a 2 R n 2 + 2 u 2 ri 

1/2 

(1 + 0.5 K)i 


1105 ( narx Y 1 

. o-Q<j rx + a 2 R n 2 + 2crA 


(1 + K)h 


\ CTrx ) 1 + 0.5 K 


(48) 
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V. Conclusions 


This paper was motivated by the desire to 1) determine what metrics should be used to characterize the 
filter performance as “good,” and 2) explore what parameters in the Kalman filter have the most impact on 
the performance of the navigation system. Closed loop control performance, which the navigation system 
exists to support, is dependent on accurate knowledge of the semimajor axis. Analytical work showed that a 
Kalman filter will not be able to produce the requisite correlation without a process noise level significantly 
lower than is likely to be possible. Further work showed that the balance and correlation requirements are 
incompatible. An examination of the effects of velocity measurements demonstrated that such measurements 
tend to reduce the correlation magnitude, making a p — > — 1 even more difficult to achieve. Increasing the 
time step of the discrete filter implementation tended to strengthen the correlation, but because the balance 
property was not simultaneously satisfied, the SMA accuracy also degraded. 
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